A microbial diversity perspective within the Restoring Resilient Ecosystems project

Dig Deeper - Microbial Community Analysis

Introduction

Biodiversity loss, ecosystem degradation, and habitat destruction are increasingly linked to human-driven changes in land use, including urbanisation, agriculture, and the exploitation of natural resources (European Parliament, 2025; Jaureguiberry et al., 2022). In response, governments across Europe — including the EU — have introduced ambitious environmental strategies such as the EU Biodiversity Strategy for 2030 (European Parliament, 2025) and the 30x30 target (Markwick, 2023), which aims to protect 30% of land and sea by the year 2030.

Ecological restoration plays a vital role in addressing these challenges. Rather than simply returning ecosystems to a previous state, modern approaches focus on restoring ecological processes and enhancing ecosystem resilience (Hicks, 2023).

The RestREco Initiative

RestREco (Restoring Resilient Ecosystems) is a NERC-funded research project that adopts a resilience-based perspective on ecological restoration. The initiative brings together researchers from:

• Cranfield University
  
• University of Stirling
• UK Centre for Ecology & Hydrology
  
• The National Trust
  
• Forest Research

Using a natural experiment design, RestREco studies a network of 133 ecological restoration sites across England and Scotland. The project aims to identify key drivers of ecosystem development, such as:

• Time since restoration began
  
• Initial ecological conditions
  
• Proximity to existing woodland and grassland

The goal is to understand how these factors influence ecosystem complexity, function, and resilience to future pressures (RestREco, 2024).

The Dig Deeper Study

As part of the RestREco initiative, the Dig Deeper study focused on how the age of restoration, establishment type, and site management affect soil microbial communities, specifically bacteria and fungi.

To explore this, high-throughput sequencing was conducted on:

• 16S rRNA gene (for bacterial communities)
  
• ITS region (for fungal communities)

The analysis focused on three main aspects:

• Alpha and beta diversity
  
• Taxonomic composition
  
• Functional diversity

These microbial assessments complement broader ecosystem-level measurements within the RestREco project, including vegetation, invertebrates, and ecosystem functions such as litter decomposition, pollination services, and soil thermodynamic efficiency.

The following sections describe the sampling design, metadata structure, and the processing pipeline used to characterise microbial communities.

Research Questions and Hypotheses

This study investigates how grassland site age, establishment method, management practices, and soil pH influence the diversity, taxonomic composition, and functional profiles of bacterial and fungal communities during restoration. It also explores the interactions between bacteria and fungi, focusing on potential correlations between microbial taxa, functional pathways, and fungal guilds.

The following hypotheses were formulated:

1. Site age, establishment method (green hay/bush (GH), seed mix (SM), natural regeneration (NR)), management type (cutting, ploughing, sheep, cattle), and pH affect alpha and beta diversity in bacterial and fungal communities.
2. These same variables influence the taxonomic composition of soil microbial communities.
3. Age, establishment method, and pH influence the functional diversity of bacterial and fungal communities.
4. There are correlations in taxon abundance between bacterial and fungal groups.
5. There are correlations in the abundance of functional pathways and guilds between bacterial and fungal communities.

---

Sample Design and Metadata Overview

Sample Collection and Geographic Coverage

A total of 330 soil samples were collected in 66 sites of England for each marker (5 per site). Although five were removed due to being incomplete.

Sample Zone - Based on GPS coordinates

Overview of the soil sampling and sequencing strategy for each microbial marker.

Metric	16S	ITS
Microbial group	Bacteria	Fungi
Region sampled	England	England
Number of sites	66	66
Samples per site	5	5
Total samples	330	330
Average reads per sample	~65,000	~65,000
Read count range	30,000–85,000	10,000–90,000

Sampling Summary

Metadata Overview

Each sample collected was accompanied by metadata capturing key environmental and management variables. These contextual factors were essential for interpreting variation in microbial diversity.

Description of metadata variables associated with soil samples

Variable	Description
Site	Name of the sampling site
Plot number	Subdivision of each site (usually 5 plots per site)
CU Code	Unique code for each sample
Year_est	Year of establishment of the site
Age	Site age (ranging from 1 year to over 100 years )
Latitude/Longitude	GPS coordinates of the sample
Establishment	Restoration type or land management
pH	Soil pH value at the time of sampling
EC	Electrical conductivity of the soil
Cutting	Whether the site is cut (1 = Yes, 0 = No)
Cattle	Presence of cattle grazing (1 = Yes, 0 = No)
Sheep	Presence of sheep grazing (1 = Yes, 0 = No)
Plough	Whether the soil has been ploughed (1 = Yes, 0 = No)

Metadata Summary

---

Analysis Pipeline

Sample Processing and Sequencing Overview

A total of 330 soil samples were collected across 66 sites (five per site) by the RestREco team. After sieving (2 mm), all samples were frozen until DNA extraction, which was performed using the QIAGEN PowerSoil kit. Amplicon sequencing targeted two regions: the V4–V5 region of the bacterial 16S rRNA gene and the ITS1 region of the fungal ITS gene. Sequencing was carried out by Novagene using Illumina paired-end reads (2×250 bp), resulting in 660 FASTQ files per marker (forward and reverse reads per sample). Metadata accompanying each sample included site characteristics (location, year of establishment, method), management practices (cutting, ploughing, grazing), and soil physicochemical properties.

16S Amplicon Processing

Quality control of 16S reads was performed using FastQC v0.12.1 and summarised with MultiQC v1.14. A custom Bash script (S02_qc.sh) ensured correct pairing of reads and removed five incomplete samples. Sequences were then imported into QIIME2 v2022.10 using the PairedEndFastqManifestPhred33V2 format, and primers were removed with CutAdapt v4.4. Read quality was visualised using qiime demux summarize.

Denoising and read pairing were performed using the DADA2 plugin (QIIME2 v2024.2), generating ASV feature tables, representative sequences, and denoising statistics. Low-abundance ASVs (<10 reads) and singletons were filtered out. Representative sequences were aligned using MAFFT and used to build phylogenetic trees with FastTree.

The ASV table was aggregated at the site level using a median-ceiling method, and rarefaction analysis was used to determine a sampling depth of 9473 reads/sample. Alpha (Observed Features, Shannon, Evenness, Faith’s PD) and beta diversity (UniFrac, Jaccard, Bray-Curtis) metrics were calculated using this rarefied table and the rooted tree. Group comparisons (e.g., by establishment, livestock, pH) were assessed using Kruskal-Wallis (alpha) and PERMANOVA (beta).

Taxonomy was assigned using a Naive Bayes classifier trained on Greengenes 13_8 (515F/806R). Taxonomic barplots were created from feature tables grouped by site. Tables were then collapsed at genus and family levels for ANCOM differential abundance testing. Pseudocounts were added where required.

Predicted functional profiles were generated using PICRUSt2 to infer KEGG Orthology (KO), EC numbers, and pathways. High-abundance pathways were retained (≥600,000 reads in ≥60 samples), and visualised with heatmaps across metadata variables. Functional differential abundance was assessed using ANCOM.

Functional alpha and beta diversity metrics were also computed, and significance was tested per metadata variable, using Kruskal-Wallis and PERMANOVA, respectively.

ITS Amplicon Processing

Fungal ITS1 reads were processed using ITSxpress v2.1.4, which performed adapter trimming and retained fungal ITS1 regions. Quality was assessed with FastQC and summarised using MultiQC v1.28 (F01_itsexpress.sh). Trimmed reads were imported into QIIME2 v2024.10, and denoised using DADA2 with no truncation.

ASV tables were aggregated by site and taxonomic assignment was performed using a pre-trained classifier based on the UNITE v10 database. Diversity analyses (Observed Features, Shannon, Evenness, Jaccard, Bray-Curtis) were conducted using a rarefied depth of 30,000 reads/sample, and visualised using PCoA and Emperor plots. Significance testing followed the same approach as for 16S (Kruskal-Wallis and PERMANOVA).

Taxonomic tables were collapsed at the genus level for ANCOM analyses. Ecological functional guilds were predicted using FUNGuild v1.1, based on collapsed species-level data formatted with custom Python scripts. Functional diversity metrics were then computed at the guild level, and visualised using core-metrics tools at a rarefied depth of 11,000 reads/sample. Significance was assessed using Kruskal-Wallis and PERMANOVA (F09 and F10 scripts).

16S–ITS Correlation Analyses

To explore inter-kingdom interactions, correlation analyses were conducted between bacterial and fungal taxa and between predicted bacterial pathways (PICRUSt2) and fungal functional guilds (FUNGuild). Using R v4.4.2 and the psych::corr.test function, Pearson correlations were calculated at the genus level. Results were adjusted using Holm-Bonferroni correction. Only strong and significant correlations (e.g., |r| > 0.7, adjusted p ≤ 0.05) were retained after filtering. Additional filters removed taxa/guilds with low overall abundance.

For functional correlations, similar thresholds were applied, but due to limited signal strength, analyses were interpreted conservatively and filtered results were reported alongside raw values when relevant.

Summary of Pipeline

Workflow (16S)

Workflow (ITS)

---

Data Pre-Processing

MutliQC on raw data

Bacteria (16S)

You can explore the full MultiQC report by clicking the image below:

MultiQC Plot (16S)

’)

Fungi (ITS)

You can explore the full MultiQC report by clicking the image below:

MultiQC Plot (ITS)

’)

QC after denoising

Statistics Table Summary

After denoising, quality control (QC) metrics provide essential insights into the effectiveness of the data processing steps and the overall quality of the resulting feature table. This summary table presents key statistics for each sample, including the number of input reads, filtered reads, and final feature counts after denoising. These metrics help assess sequencing success, identify potential outliers, and ensure that sufficient data remain for robust downstream analyses. Samples with unusually low read counts or feature richness may need to be excluded or interpreted with caution.

Bacteria (16S)

Here is a link to the statistics after denoising to view it on QIIME2 (16S) : Statitics after denoising (16S)

Statitics after denoising (16S)

Fungi (ITS)

Here is a link to the statistics after denoising to view it on QIIME2 (ITS) : Statitics after denoising (16S)

Statitics after denoising (ITS)

QC Plots

These Quality Control (QC) plots were generated after trimming the sequencing reads. They provide a visual summary of the base quality scores, read length distributions, and other metrics, helping to assess whether the trimming step successfully removed low-quality regions and adapter contamination. Consistently high-quality reads across samples are essential for reliable downstream analysis.

Bacteria (16S)

QC plot (16S)

Fungi (ITS)

QC plot (ITS)

Rarefication Curves

Rarefaction curves provide a visual tool to assess sequencing depth and compare species richness between samples. In this study, rarefaction curves were generated separately for bacterial (16S rRNA gene) and fungal (ITS) communities, using the number of observed features—i.e., unique ASVs—as a proxy for richness.

For both bacteria and fungi, the shape of each curve indicates whether sequencing depth was sufficient to capture most of the diversity in a given sample. Curves that level off suggest that a representative portion of the community has been sampled, whereas rising curves indicate that additional sequencing could reveal further diversity.

This step is crucial to ensure that downstream diversity analyses are not biased by unequal sampling effort.

Bacteria (16S)

Here is a link to the Rarefiction plots for more flexibility on QIIME2: Rarefiction plots (16S)

Select variable: Establishment Age Year of Establishment Os Location Plough Lat/long

Rarefaction Curves of Observed Features by 'Establishment' (16S)

Rarefaction Curves of Observed Features by 'Age' (16S)

Rarefaction Curves of Observed Features by 'Year_est' (16S)

Rarefaction Curves of Observed Features by 'OS_location' (16S)

Rarefaction Curves of Observed Features by 'Lat_long' (16S)

Rarefaction Curves of Observed Features by 'Plough' (16S)

---

Fungi (ITS)

Here is a link to the Rarefiction plots for more flexibility on QIIME2: Rarefiction plots (ITS)

Select variable: Os Location CU Code Lat/long

Rarefaction Curves of Observed Features by 'OS_location' (ITS)

Rarefaction Curves of Observed Features by 'CU Code' (ITS)

Rarefaction Curves of Observed Features by 'Lat_long' (ITS)

---

Taxonomic Diversity (Bacteria & Fungi)

Results Summary
Microbial taxonomic diversity revealed clear responses to restoration strategies. For bacteria, establishment method was the strongest driver of alpha and beta diversity, followed by pH and grazing. For fungi, pH and age were more influential. NR sites generally showed lower bacterial diversity and functional richness than GH and SM.

Alpha Diversity

Alpha diversity refers to the variety of organisms within a particular sample or environment. It reflects both richness—the number of distinct taxa—and evenness—how evenly individuals are distributed among those taxa.

To allow interactive exploration of alpha diversity metrics across different environmental variables, we implemented a drop-down menu that dynamically displays the corresponding plots. Some variables, such as pH category, are only present in the ITS dataset, while others, like Year group, are specific to the 16S dataset. Internally, variables are mapped to their dataset-specific equivalents where needed (e.g. Age group in 16S becomes Age category in ITS). It is important to note, however, that these variables are not always directly comparable: for instance, Age group (16S) divides sites into multiple discrete intervals based on restoration age, while Age category (ITS) is a binary classification based on whether a site is above or below the median age ( Age_binary for bacteria). Despite these differences, the interface ensures that only available and relevant plots are shown for each selection.

Bacteria (16S)

Results Summary

Shannon entropy revealed significant overall differences, with green hay/bush sites showing higher diversity than both natural regeneration and seed mix sites. Seed mix sites also had higher Shannon diversity than natural regeneration. Pielou’s evenness differed significantly across establishment methods (p = 0.0148), with green hay/bush sites being significantly more even than natural regeneration (q = 0.0148). The difference between seed mix and natural regeneration approached significance (q = 0.0597), while no significant difference was found between green hay/bush and seed mix.

Select variable: Establishment Cutting Cattle Sheep Plough Year group Age group pH binary Age binary

Shannon Index Boxplots

Click for more information on Shannon index

The Shannon index takes into account not only the number of species present, but also how evenly their abundances are distributed. A higher Shannon value generally indicates a more diverse and ecologically balanced community.

Here is a link to the full QIIME2 results (16S) : Shannon Index (16S)

Kruskal-Wallis p-value: 0.00102

Kruskal-Wallis p-value: 0.158

Kruskal-Wallis p-value: 0.939

Kruskal-Wallis p-value: 0.655

Kruskal-Wallis p-value: 0.316

Kruskal-Wallis p-value: 0.408

Kruskal-Wallis p-value: 0.126

Kruskal-Wallis p-value: 0.116

Kruskal-Wallis p-value: 0.358

Faith PD Boxplots

Click for more information on Faith PD

Faith’s Phylogenetic Diversity (Faith PD) measures the total branch length of the phylogenetic tree that spans the species in a sample. Unlike the Shannon index, Faith PD incorporates evolutionary relationships, providing a phylogenetic perspective on diversity. (Bacteria only)

Here is a link to the full QIIME2 results (16S) : Faith PD (16S)

Kruskal-Wallis p-value: 0.0106

Kruskal-Wallis p-value: 0.362

Kruskal-Wallis p-value: 0.878

Kruskal-Wallis p-value: 0.604

Kruskal-Wallis p-value: 0.705

Kruskal-Wallis p-value: 0.216

Kruskal-Wallis p-value: 0.0709

Kruskal-Wallis p-value: 0.0893

Kruskal-Wallis p-value: 0.393

Evenness Index Boxplots

Click for more information on Pielou’s Evenness Index

Pielou’s Evenness index specifically quantifies how equally individual organisms are distributed across taxa. While Shannon integrates both richness and evenness, this metric isolates the evenness component, providing a complementary view of diversity patterns.

Here is a link to the full QIIME2 results (16S) : Pielou Evenness (16S)

Kruskal-Wallis p-value: 0.0148

Kruskal-Wallis p-value: 0.17

Kruskal-Wallis p-value: 0.878

Kruskal-Wallis p-value: 0.291

Kruskal-Wallis p-value: 0.64

Kruskal-Wallis p-value: 0.854

Kruskal-Wallis p-value: 0.608

Kruskal-Wallis p-value: 0.542

Kruskal-Wallis p-value: 0.0817

---

Fungi (ITS)

For the variable pH_category, the median is 7.94 and for Age_category it’s 13 years.

Results Summary

Fungal alpha diversity showed limited variation across establishment methods, with natural regeneration sites displaying slightly higher diversity, though differences were not statistically significant. Similarly, grazing, ploughing, and cutting had no detectable impact on overall fungal diversity.

However, soil pH significantly affected fungal diversity, with alkaline soils (pH > 7.9) showing lower Shannon entropy values compared to neutral or slightly acidic soils (p < 0.02). Site age also influenced diversity, with older sites displaying lower median Shannon diversity (p < 0.02).

Regarding community evenness, sheep grazing was the only factor associated with a significant difference: grazed sites had slightly higher Pielou’s evenness (p = 0.03). Other variables, including establishment method, age, and pH, showed no significant effect on evenness.

Select variable: Establishment Cutting Cattle Sheep Plough Age category pH category

Shannon Index Boxplots

Click for more information on Shannon index

The Shannon index takes into account not only the number of species present, but also how evenly their abundances are distributed. A higher Shannon value generally indicates a more diverse and ecologically balanced community.

Here is a link to the full QIIME2 results (ITS) : Shannon Index (ITS)

Kruskal-Wallis p-value: 0.796

Kruskal-Wallis p-value: 0.719

Kruskal-Wallis p-value: 0.768

Kruskal-Wallis p-value: 0.308

Kruskal-Wallis p-value: 0.647

Kruskal-Wallis p-value: 0.159

Kruskal-Wallis p-value: 0.193

Evenness Index Boxplots

Click for more information on Pielou’s Evenness Index

Pielou’s Evenness index specifically quantifies how equally individual organisms are distributed across taxa. While Shannon integrates both richness and evenness, this metric isolates the evenness component, providing a complementary view of diversity patterns.

Here is a link to the full QIIME2 results (ITS) : Pielou Evenness (ITS)

Kruskal-Wallis p-value: 0.512

Kruskal-Wallis p-value: 0.719

Kruskal-Wallis p-value: 0.928

Kruskal-Wallis p-value: 0.0927

Kruskal-Wallis p-value: 0.389

Kruskal-Wallis p-value: 0.386

Kruskal-Wallis p-value: 0.954

---

Comparative Microbial Community Composition (Beta Diversity)

To explore differences in microbial communities, we often rely on dimensionality reduction techniques such as Principal Coordinates Analysis (PCoA), visualised through Emperor plots. Two commonly used distance metrics in this context are Bray-Curtis and Jaccard.

While both metrics can reveal meaningful clustering and separation in microbial data, they capture complementary aspects of community structure.

This approach is useful to visualise group clustering by variables like establishment type, management, or site age, and complements statistical tests like PERMANOVA, which assess whether community structure significantly varies across those factors.

Results Summary

Beta diversity patterns confirmed that establishment type, pH, and grazing significantly shaped microbial community composition. For bacteria, strong clustering by method and significant PERMANOVA results support distinct assemblages. Fungal beta diversity followed similar trends but with additional sensitivity to site age.

Bacteria (16S)

Bray-Curtis

Click for more information on Bray-Curtis Distance Metric

Bray-Curtis is a non-phylogenetic beta diversity metric that compares microbial communities based on species abundance, taking into account both shared taxa and their relative proportions. Unlike presence/absence metrics, Bray-Curtis gives more weight to dominant species and is sensitive to changes in abundance structure. It ranges from 0 (identical composition and abundance) to 1 (completely different), reflecting how much two communities differ in the quantity of each taxon they contain. Bray-Curtis is especially useful for detecting shifts in community structure driven by changes in resource availability, disturbance, or other factors that influence the balance of taxa rather than just their presence.

Results Summary

The Bray–Curtis metric revealed strong effects of establishment method, soil pH, and sheep grazing on bacterial community structure. Significant differences were detected across all establishment types (PERMANOVA p = 0.001), with clear clustering observed in PCoA space. Soil age also influenced community structure, with younger and older sites differing significantly (p = 0.035), despite weak visual separation.

Emperor Plot

Results Summary

Beta diversity analysis using Bray–Curtis distance revealed that establishment method, soil pH, and sheep grazing were the primary factors shaping bacterial community structure. Significant differences in taxonomic composition were observed between all establishment types, with natural regeneration, seed mix, and green hay/bush plots forming distinct but partially overlapping clusters in PCoA space. Site age also had a detectable effect with Bray–Curtis, highlighting changes in relative abundance between younger and older sites.

Here is a link to the Bray-Curtis Emperor Plot for more flexibility on QIIME2: Bray-Curtis Emperor Plot (16S)

Bray-Curtis Emperor Plot

PERMANOVA

Results Summary

Bray-Curtis-based PERMANOVA revealed significant differences in bacterial community composition across establishment methods, with all pairwise comparisons between methods also significant. These results indicate that restoration strategy has a strong influence on both the composition and abundance structure of bacterial communities.

Site age also had a significant but weaker effect. However, clustering by age in Bray-Curtis PCoA plots was less pronounced, suggesting more subtle shifts in community structure over time.

Soil pH significantly shaped bacterial community composition, with Bray-Curtis PCoA plots showing clear clustering between high and low-to-neutral pH soils

You can click on the images below to access the full QIIME2 report.

Select variable: Establishment Cutting Cattle Sheep Plough Year group pH binary Age binary

View full QIIME2 results (Bray–Curtis – Establishment)

Figure: PERMANOVA Bray–Curtis for “Establishment” (16S)

View full QIIME2 results (Bray–Curtis – Cutting)

Figure: PERMANOVA Bray–Curtis for “Cutting” (16S)

View full QIIME2 results (Bray–Curtis – Cattle)

Figure: PERMANOVA Bray–Curtis for “Cattle” (16S)

View full QIIME2 results (Bray–Curtis – Sheep)

Figure: PERMANOVA Bray–Curtis for “Sheep” (16S)

View full QIIME2 results (Bray–Curtis – Plough)

Figure: PERMANOVA Bray–Curtis for “Plough” (16S)

View full QIIME2 results (Bray–Curtis – Year_group)

Figure: PERMANOVA Bray–Curtis for “Year_group” (16S)

View full QIIME2 results (Bray–Curtis – pH_binary)

Figure: PERMANOVA Bray–Curtis for “pH_binary” (16S)

View full QIIME2 results (Bray–Curtis – Age_binary)

Figure: PERMANOVA Bray–Curtis for “Age_binary” (16S)

Jaccard

Click for more information on Jaccard Distance Metric

Jaccard is a non-phylogenetic beta diversity metric that compares microbial communities based on the presence or absence of taxa, without considering their abundance or evolutionary relationships. It measures the proportion of shared taxa between samples relative to the total number of taxa present, treating all species equally. This makes it particularly sensitive to community turnover, highlighting whether two sites contain the same organisms, regardless of how abundant they are. Jaccard is well-suited for detecting changes in community membership across environmental gradients or treatment groups, especially when interested in whether taxa are present at all, rather than how dominant they are.

Emperor Plot

Here is a link to the Jaccard Emperor Plot for more flexibility on QIIME2: Jaccard Emperor Plot (16S)

Jaccard Emperor Plot (16S)

PERMANOVA

You can click on the images below to access the full QIIME2 report.

Select variable: Establishment Cutting Cattle Sheep Plough Year_group pH binary Age binary

View full QIIME2 results (Jaccard – Establishment)

Figure: PERMANOVA Jaccard for “Establishment” (16S)

View full QIIME2 results (Jaccard – Cutting)

Figure: PERMANOVA Jaccard for “Cutting” (16S)

View full QIIME2 results (Jaccard – Cattle)

Figure: PERMANOVA Jaccard for “Cattle” (16S)

View full QIIME2 results (Jaccard – Sheep)

Figure: PERMANOVA Jaccard for “Sheep” (16S)

View full QIIME2 results (Jaccard – Plough)

Figure: PERMANOVA Jaccard for “Plough” (16S)

View full QIIME2 results (Jaccard – Year_group)

Figure: PERMANOVA Jaccard for “Year_group” (16S)

View full QIIME2 results (Jaccard – pH_binary)

Figure: PERMANOVA Jaccard for “pH_binary” (16S)

View full QIIME2 results (Jaccard – Age_binary)

Figure: PERMANOVA Jaccard for “Age_binary” (16S)

Unweighted Unifrac

Click for more information on Unweighted Unifrac Distance Metric

Unweighted UniFrac is a phylogenetic beta diversity metric that compares microbial communities based solely on the presence or absence of taxa, while incorporating their evolutionary relationships.

Unlike abundance-based measures, this metric considers whether lineages are shared between communities, regardless of how dominant they are. It is particularly sensitive to rare or low-abundance taxa, as all taxa are weighted equally.

Unweighted UniFrac is useful for identifying broad differences in community membership—for example, whether two sites share the same species, even if those species occur at very different abundances. This makes it well-suited to detecting compositional shifts due to strong environmental filters or historical legacies.

Results Summary

Unweighted UniFrac, which focuses on phylogenetic presence/absence, similarly highlighted establishment method, soil pH, and sheep grazing as major drivers of variation. All pairwise comparisons between establishment types yielded statistically significant differences (q < 0.005), confirming shifts in community membership across restoration strategies. However, site age had no significant effect using this metric (p = 0.063), suggesting that observed changes with age were driven more by shifts in abundance than by changes in phylogenetic composition.

Emperor Plot

Here is a link to the Unweighted Unifrac Emperor Plot for more flexibility on QIIME2: Unweighted Unifrac Emperor Plot (16S)

Unweighted Unifrac Emperor Plot (16S)

PERMANOVA

You can click on the images below to access the full QIIME2 report.

Select variable: Establishment Cutting Cattle Sheep Plough Year_group pH binary Age binary

View full QIIME2 results (Unweighted Unifrac – Establishment)

Figure: PERMANOVA Unweighted Unifrac for “Establishment” (16S)

View full QIIME2 results (Unweighted Unifrac – Cutting)

Figure: PERMANOVA Unweighted Unifrac for “Cutting” (16S)

View full QIIME2 results (Unweighted Unifrac – Cattle)

Figure: PERMANOVA Unweighted Unifrac for “Cattle” (16S)

View full QIIME2 results (Unweighted Unifrac – Sheep)

Figure: PERMANOVA Unweighted Unifrac for “Sheep” (16S)

View full QIIME2 results (Unweighted Unifrac – Plough)

Figure: PERMANOVA Unweighted Unifrac for “Plough” (16S)

View full QIIME2 results (Unweighted Unifrac – Year_group)

Figure: PERMANOVA Unweighted Unifrac for “Year_group” (16S)

View full QIIME2 results (Unweighted Unifrac – pH_binary)

Figure: PERMANOVA Unweighted Unifrac for “pH_binary” (16S)

View full QIIME2 results (Unweighted Unifrac – Age_binary)

Figure: PERMANOVA Unweighted Unifrac for “Age_binary” (16S)

Weighted Unifrac

Click for more information on Weighted Unifrac Distance Metric

Weighted UniFrac is a phylogenetic distance metric that takes into account both the evolutionary relationships between taxa and their relative abundances in each sample.

It quantifies how much of the phylogenetic tree is shared between communities, with branches weighted by the proportion of abundance they represent. As a result, changes in dominant taxa have more influence on the distance than rare species.

Weighted UniFrac is ideal for detecting differences in community structure, especially when those differences involve shifts in abundant lineages. It provides a more nuanced view than unweighted metrics by integrating both phylogenetic and quantitative information.

Emperor Plot

Here is a link to the Weighted Unifrac Emperor Plot for more flexibility on QIIME2: Weighted Unifrac Emperor Plot (16S)

Weighted Unifrac Emperor Plot (16S)

PERMANOVA

You can click on the images below to access the full QIIME2 report.

Select variable: Establishment Cutting Cattle Sheep Plough Year_group pH binary Age binary

View full QIIME2 results (Weighted Unifrac – Establishment)

Figure: PERMANOVA Weighted Unifrac for “Establishment” (16S)

View full QIIME2 results (Weighted Unifrac – Cutting)

Figure: PERMANOVA Weighted Unifrac for “Cutting” (16S)

View full QIIME2 results (Weighted Unifrac – Cattle)

Figure: PERMANOVA Weighted Unifrac for “Cattle” (16S)

View full QIIME2 results (Weighted Unifrac – Sheep)

Figure: PERMANOVA Weighted Unifrac for “Sheep” (16S)

View full QIIME2 results (Weighted Unifrac – Plough)

Figure: PERMANOVA Weighted Unifrac for “Plough” (16S)

View full QIIME2 results (Weighted Unifrac – Year_group)

Figure: PERMANOVA Weighted Unifrac for “Year_group” (16S)

View full QIIME2 results (Weighted Unifrac – pH_binary)

Figure: PERMANOVA Weighted Unifrac for “pH_binary” (16S)

View full QIIME2 results (Weighted Unifrac – Age_binary)

Figure: PERMANOVA Weighted Unifrac for “Age_binary” (16S)

---

Fungi (ITS)

Bray-Curtis

Click for more information on Bray-Curtis Distance Metric

Bray-Curtis is a non-phylogenetic beta diversity metric that compares microbial communities based on species abundance, taking into account both shared taxa and their relative proportions. Unlike presence/absence metrics, Bray-Curtis gives more weight to dominant species and is sensitive to changes in abundance structure. It ranges from 0 (identical composition and abundance) to 1 (completely different), reflecting how much two communities differ in the quantity of each taxon they contain. Bray-Curtis is especially useful for detecting shifts in community structure driven by changes in resource availability, disturbance, or other factors that influence the balance of taxa rather than just their presence.

Emperor Plot

Results Summary

Beta diversity analysis revealed that establishment method, soil pH, and site age had the most significant influence on fungal community composition, followed by sheep grazing, ploughing, and cutting. In PCoA plots, natural regeneration sites formed a loose cluster near seed mix sites, while green hay/bush sites clustered more distinctly along Axis 1.

Here is a link to the Bray-Curtis Emperor Plot for more flexibility on QIIME2: Bray-Curtis Emperor Plot (ITS)

Bray-Curtis Emperor Plot (ITS)

PERMANOVA

Results Summary

PERMANOVA analysis identified pH category (p = 0.001) and establishment method (p = 0.002) as the strongest drivers of differences in fungal communities. Site age (p = 0.024) and sheep grazing (p = 0.025) also had significant effects.

Pairwise comparisons showed the strongest differences between green hay/bush and seed mix sites (p = 0.008), followed by green hay/bush vs. natural regeneration (p = 0.011). Natural regeneration sites exhibited the highest variation in community composition, with Bray-Curtis distances ranging from 0.25 to 0.95.

You can click on the images below to access the full QIIME2 report.

Select variable: Establishment Cutting Cattle Sheep Plough pH category Age category

View full QIIME2 results (Bray–Curtis – Establishment)

Figure: PERMANOVA Bray–Curtis for “Establishment” (ITS)

View full QIIME2 results (Bray–Curtis – Cutting)

Figure: PERMANOVA Bray–Curtis for “Cutting” (ITS)

View full QIIME2 results (Bray–Curtis – Cattle)

Figure: PERMANOVA Bray–Curtis for “Cattle” (ITS)

View full QIIME2 results (Bray–Curtis – Sheep)

Figure: PERMANOVA Bray–Curtis for “Sheep” (ITS)

View full QIIME2 results (Bray–Curtis – Plough)

Figure: PERMANOVA Bray–Curtis for “Plough” (ITS)

View full QIIME2 results (Bray–Curtis – pH_category)

Figure: PERMANOVA Bray–Curtis for “pH_category” (ITS)

View full QIIME2 results (Bray–Curtis – Age_category)

Figure: PERMANOVA Bray–Curtis for “Age_category” (ITS)

Jaccard

Click for more information on Jaccard Distance Metric

Jaccard is a non-phylogenetic beta diversity metric that compares microbial communities based on the presence or absence of taxa, without considering their abundance or evolutionary relationships. It measures the proportion of shared taxa between samples relative to the total number of taxa present, treating all species equally. This makes it particularly sensitive to community turnover, highlighting whether two sites contain the same organisms, regardless of how abundant they are. Jaccard is well-suited for detecting changes in community membership across environmental gradients or treatment groups, especially when interested in whether taxa are present at all, rather than how dominant they are.

Emperor Plot

Results Summary

Jaccard analysis confirmed similar drivers of community structure—establishment method, pH, and age—with natural regeneration, seed mix, and green hay/bush sites forming distinguishable but overlapping clusters in higher dimensions. This suggests that both presence/absence and abundance-based differences in fungal communities are shaped by restoration strategy and site conditions.

Here is a link to the Jaccard Emperor Plot for more flexibility on QIIME2: Jaccard Emperor Plot (ITS)

Jaccard Emperor Plot (ITS)

PERMANOVA

Results Summary

With Jaccard dissimilarity, establishment method (p = 0.001) and pH (p = 0.002) remained the strongest drivers, followed by site age (p = 0.006), ploughing (p = 0.024), and cutting (p = 0.04).

Pairwise comparisons showed the most significant differences between green hay/bush and both natural regeneration and seed mix sites (p = 0.001), and also between natural regeneration and seed mix (p = 0.004). Median Jaccard dissimilarity values were high across all establishment methods (0.75–0.85), indicating substantial variation in taxon presence/absence.

You can click on the images below to access the full QIIME2 report.

Select variable: Establishment Cutting Cattle Sheep Plough pH category Age category

View full QIIME2 results (Jaccard – Establishment)

Figure: PERMANOVA Jaccard for “Establishment” (ITS)

View full QIIME2 results (Jaccard – Cutting)

Figure: PERMANOVA Jaccard for “Cutting” (ITS)

View full QIIME2 results (Jaccard – Cattle)

Figure: PERMANOVA Jaccard for “Cattle” (ITS)

View full QIIME2 results (Jaccard – Sheep)

Figure: PERMANOVA Jaccard for “Sheep” (ITS)

View full QIIME2 results (Jaccard – Plough)

Figure: PERMANOVA Jaccard for “Plough” (ITS)

View full QIIME2 results (Jaccard – pH_category)

Figure: PERMANOVA Jaccard for “pH_category” (ITS)

View full QIIME2 results (Jaccard – Age_category)

Figure: PERMANOVA Jaccard for “Age_category” (ITS)

---

Taxonomy Composition

Taxonomy Barplot

Bacteria (16S)

Results Summary

Analysis of bacterial community composition at the phylum level across establishment types revealed that Actinobacteria and Proteobacteria were dominant, together accounting for 50–80% of relative abundance. Actinobacteria were particularly abundant in seed mix (SM) sites, where they reached up to ~48%, while their proportion was lowest in natural regeneration (NR) sites (~21%). In contrast, Proteobacteria were most prominent in NR soils, peaking at ~55%. Other frequently detected phyla included Acidobacteria, Firmicutes, Chloroflexi, Verrucomicrobia, Planctomycetes, and Nitrospirae, with Verrucomicrobia notably more abundant in NR plots.

At the class level, Alphaproteobacteria, Actinobacteria, Thermoleophilia, and Bacilli dominated across sites. A higher proportion of Spartobacteria was observed in NR plots compared to GH and SM plots.

When comparing soils with pH above and below the median, the dominant phyla remained broadly consistent. However, some shifts were noted: Firmicutes and Planctomycetes were more abundant in more acidic soils (below-median pH), with Firmicutes reaching up to ~30% (versus ~19% above the median), and Planctomycetes peaking at ~6.5% (compared to ~3.6%).

Here is a link to the Taxonomy Barplots for more flexibility on QIIME2: Taxonomy Barplot (16S)

Select variable: Establishment Cutting Cattle Sheep Plough Age group pH binary Year group Age binary

Taxonomy Barplots Associated with 'Establishment' (16S)

Taxonomy Barplots Associated with 'Cutting' (16S)

Taxonomy Barplots Associated with 'Cattle' (16S)

Taxonomy Barplots Associated with 'Sheep' (16S)

Taxonomy Barplots Associated with 'Plough' (16S)

Taxonomy Barplots Associated with 'Age_group' (16S)

Taxonomy Barplots Associated with 'pH_binary' (16S)

Taxonomy Barplots Associated with 'Year_group' (16S)

Taxonomy Barplots Associated with 'Age_binary' (16S)

Fungi (ITS)

Results Summary

Across most sites, the dominant fungal phyla were Ascomycota, Mortierellomycota, and Basidiomycota, alongside a notable proportion of unclassified fungi (incertae sedis). Minor differences in community composition were observed among natural regeneration (NR) plots, with Basidiomycota being more prominent in one site.

At the class level, Sordariomycetes, Mortierellomycetes, Dothideomycetes, and Leotiomycetes were commonly detected across all sites. However, Agaricomycetes appeared less frequently in NR plots compared to GH and SM sites, except in one case.

When grouping sites by pH category, the most visible difference was an increased presence of Sordariomycetes in slightly acidic to neutral soils (below the median pH), compared to more alkaline environments.

Here is a link to the Taxonomy Barplots for more flexibility on QIIME2: Taxonomy Barplot (ITS)

Select variable: Establishment Cutting Cattle Sheep Plough Age category pH category

Taxonomy Barplots Associated with 'Establishment' (ITS)

Taxonomy Barplots Associated with 'Cutting' (ITS)

Taxonomy Barplots Associated with 'Cattle' (ITS)

Taxonomy Barplots Associated with 'Sheep' (ITS)

Taxonomy Barplots Associated with 'Plough' (ITS)

Taxonomy Barplots Associated with 'Age_category' (ITS)

Taxonomy Barplots Associated with 'pH_category' (ITS)

Krona Plots

To explore the composition of soil microbial communities, we used Krona plots — interactive, circular charts that display taxonomic abundances in a hierarchical manner.

These plots allow users to intuitively navigate from broader taxonomic levels (such as Phylum) to more specific ones (like Genus), while simultaneously comparing relative abundances across taxa.

In this study, Krona plots provide a powerful and user-friendly way to:

• Visualise which microbial groups dominate each site
• Explore the taxonomic diversity present in bacterial, archaeal, and fungal communities

You can click on the images below to access the Krona plots for each site.

Bacteria (16S)

Krona Plot for Baltic_farm_1 (16S)

Fungi (ITS)

Krona Plot for Baltic_farm_1 (ITS)

Differential Abundance Analysis with ANCOM

We used ANCOM to identify taxa whose relative abundances significantly differed across groups. This method accounts for the compositional nature of microbiome data by comparing log-ratios between taxa. The results are shown as volcano-like plots, where the W statistic reflects how many pairwise comparisons a taxon was found to differ in. Significant taxa are highlighted accordingly.

Bacteria (16S)

Results Summary

ANCOM analysis identified several bacterial taxa whose abundance varied significantly across environmental factors.

Conexibacteraceae was more abundant in soils with below-median pH, indicating a preference for more acidic conditions.

The genus Tetrasphaera showed significantly higher abundance in younger sites, suggesting a potential role in early stages of soil restoration.

Two taxa were strongly associated with establishment method: an unclassified genus within the Intrasporangiaceae family and the genus Blastococcus were both more abundant in SM and GH plots compared to NR. These differences were supported by large W statistics and differences in median abundance values across treatments. These taxa responded notably to restoration strategy, indicating that establishment method can drive specific shifts in bacterial community structure.

Select variable: Establishment Cutting Cattle Sheep Plough Age_group pH_binary Year_group Age_binary

View interactive volcano-like plot in QIIME2 (16S - Establishment )

Differences in Taxa Abundance Associated with 'Establishment' (ANCOM Results) (16S)

View interactive volcano-like plot in QIIME2 (16S - Cutting )

Differences in Taxa Abundance Associated with 'Cutting' (ANCOM Results) (16S)

View interactive volcano-like plot in QIIME2 (16S - Cattle )

Differences in Taxa Abundance Associated with 'Cattle' (ANCOM Results) (16S)

View interactive volcano-like plot in QIIME2 (16S - Sheep )

Differences in Taxa Abundance Associated with 'Sheep' (ANCOM Results) (16S)

View interactive volcano-like plot in QIIME2 (16S - Plough )

Differences in Taxa Abundance Associated with 'Plough' (ANCOM Results) (16S)

View interactive volcano-like plot in QIIME2 (16S - Age_group )

Differences in Taxa Abundance Associated with 'Age_group' (ANCOM Results) (16S)

View interactive volcano-like plot in QIIME2 (16S - pH_binary )

Differences in Taxa Abundance Associated with 'pH_binary' (ANCOM Results) (16S)

View interactive volcano-like plot in QIIME2 (16S - Year_group )

Differences in Taxa Abundance Associated with 'Year_group' (ANCOM Results) (16S)

View interactive volcano-like plot in QIIME2 (16S - Age_binary )

Differences in Taxa Abundance Associated with 'Age_binary' (ANCOM Results) (16S)

Fungi (ITS)

Results Summary

ANCOM analysis identified several fungal taxa as differentially abundant in relation to soil pH, site age, and establishment method:

• Metapochonia (W = 386) was significantly more abundant in acidic soils, suggesting a preference for low-pH environments.

• Two taxa were age-responsive: a Lasiosphaeriaceae taxon (W = 366) was more abundant in younger soils, while Paraphaeosphaeria (W = 320) was enriched in older soils, indicating shifts in fungal community structure with site maturation.

• Regarding establishment method, Gibellulopsis (W = 596) was predominantly found in green hay/bush and seed mix sites, but rare in natural regeneration sites. In contrast, Metarhizium (W = 535) was most abundant in naturally regenerating sites, likely reflecting adaptation to less-disturbed, self-organising ecosystems.

These findings highlight how restoration strategy, pH, and time since restoration can drive specific taxonomic responses, shaping fungal community composition in grassland soils.

Select variable: Establishment Cutting Cattle Sheep Plough Age_category pH_category

View interactive volcano-like plot in QIIME2 (ITS - Establishment )

Differences in Taxa Abundance Associated with 'Establishment' (ANCOM Results) (ITS)

View interactive volcano-like plot in QIIME2 (ITS - Cutting )

Differences in Taxa Abundance Associated with 'Cutting' (ANCOM Results) (ITS)

View interactive volcano-like plot in QIIME2 (ITS - Cattle )

Differences in Taxa Abundance Associated with 'Cattle' (ANCOM Results) (ITS)

View interactive volcano-like plot in QIIME2 (ITS - Sheep )

Differences in Taxa Abundance Associated with 'Sheep' (ANCOM Results) (ITS)

View interactive volcano-like plot in QIIME2 (ITS - Plough )

Differences in Taxa Abundance Associated with 'Plough' (ANCOM Results) (ITS)

View interactive volcano-like plot in QIIME2 (ITS - Age_category )

Differences in Taxa Abundance Associated with 'Age_category' (ANCOM Results) (ITS)

View interactive volcano-like plot in QIIME2 (ITS - pH_category )

Differences in Taxa Abundance Associated with 'pH_category' (ANCOM Results) (ITS)

---

Functional Diversity (Bacteria & Fungi)

Functional diversity refers to the variety of biological functions present in microbial communities. Unlike taxonomic analysis, which identifies who is there, functional analysis explores what they can do. This section includes two parallel analyses:

• Bacteria (16S): Functional profiles were inferred using PICRUSt2, predicting gene family and pathway abundance.
• Fungi (ITS): Functional roles were explored using FUNGuild, which assigns taxa to ecological guilds based on literature.

Alpha and beta diversity were analysed, along with differentially abundant functions or pathways.

Results Summary

Functional diversity revealed strong effects of establishment method and pH on microbial pathway composition. Core metabolic functions were consistent across sites, but specific biosynthetic and degradation pathways differed by age, pH, and grazing. Saprotrophic functions were dominant in fungi, with guild shifts linked to disturbance and soil chemistry.

Bacteria (16S)

Results Summary

Bacterial functional alpha and beta diversity varied significantly across establishment types, with GH and SM plots showing enhanced diversity. Functional profiles also responded to pH and age. Several metabolic pathways, including heme and amino acid biosynthesis, were enriched in specific conditions, reflecting adaptive responses to restoration.

Functional Alpha Diversity

We examined three complementary metrics: Shannon diversity (richness and evenness), Observed Features (richness only), and Evenness (dominance balance). You can use the menu below to explore how each metric varies according to different experimental variables.

Results Summary

Functional alpha diversity, based on PICRUSt2-inferred metabolic pathways, varied significantly across establishment types. Diversity was higher in SM and GH plots compared to NR, while no significant difference was detected between SM and GH.

Select variable: Establishment Cutting Cattle Sheep Plough Age_group pH_binary Age_binary

Shannon Index Boxplots

Here is a link to the Shannon Index Boxplot for more flexibility on QIIME2: Shannon Index Boxplot - Functional (16S)

Kruskal-Wallis p-value: 0.00164

Kruskal-Wallis p-value: 0.0496

Kruskal-Wallis p-value: 0.739

Kruskal-Wallis p-value: 0.681

Kruskal-Wallis p-value: 0.588

Kruskal-Wallis p-value: 0.144

Kruskal-Wallis p-value: 0.517

Kruskal-Wallis p-value: 0.0887

Observed Features Boxplots

Here is a link to the Observed features Boxplot for more flexibility on QIIME2: Observed Features Boxplot - Functional (16S)

Kruskal-Wallis p-value: 0.0354

Kruskal-Wallis p-value: 0.758

Kruskal-Wallis p-value: 0.969

Kruskal-Wallis p-value: 0.707

Kruskal-Wallis p-value: 0.923

Kruskal-Wallis p-value: 0.372

Kruskal-Wallis p-value: 0.292

Kruskal-Wallis p-value: 0.974

Evenness Index Boxplots

Here is a link to the Evenness Boxplot for more flexibility on QIIME2: Evenness Boxplot - Functional (16S)

Kruskal-Wallis p-value: 0.0294

Kruskal-Wallis p-value: 0.13

Kruskal-Wallis p-value: 0.888

Kruskal-Wallis p-value: 0.858

Kruskal-Wallis p-value: 0.393

Kruskal-Wallis p-value: 0.84

Kruskal-Wallis p-value: 0.184

Kruskal-Wallis p-value: 0.235

Functional Beta Diversity

Emperor Plot – Bray-Curtis Distance

Results Summary

Functional beta diversity analyses based on PICRUSt2 pathway profiles showed that multiple site-level factors significantly influenced microbial community composition. Establishment method had the strongest effect, with distinct functional profiles observed across all establishment types. Soil pH also shaped functional composition, with samples from high- and low-pH sites showing partial separation. In addition, grazing by sheep was associated with significant shifts in predicted functional pathways, despite some overlap between grazed and ungrazed plots.

Here is a link to the Bray-Curtis Emperor Plot for more flexibility on QIIME2: Bray-Curtis Emperor Plot (16S)

Bray-Curtis Emperor Plot (Functional)

PERMANOVA – Bray-Curtis Distance

You can click on the images below to access the full QIIME2 report.

Select variable: Establishment Cutting Cattle Sheep Plough Year group Age group Age Binary pH Binary

View full QIIME2 results (Bray–Curtis – Establishment)

Figure: PERMANOVA Bray-Curtis for “Establishment” (ITS)

View full QIIME2 results (Bray–Curtis – Cutting)

Figure: PERMANOVA Bray-Curtis for “Cutting” (ITS)

View full QIIME2 results (Bray–Curtis – Cattle)

Figure: PERMANOVA Bray-Curtis for “Cattle” (ITS)

View full QIIME2 results (Bray–Curtis – Sheep)

Figure: PERMANOVA Bray-Curtis for “Sheep” (ITS)

View full QIIME2 results (Bray–Curtis – Plough)

Figure: PERMANOVA Bray-Curtis for “Plough” (ITS)

View full QIIME2 results (Bray–Curtis – Year_group)

Figure: PERMANOVA Bray-Curtis for “Year_group” (ITS)

View full QIIME2 results (Bray–Curtis – Age_group)

Figure: PERMANOVA Bray-Curtis for “Age_group” (ITS)

View full QIIME2 results (Bray–Curtis – Age_binary)

Figure: PERMANOVA Bray-Curtis for “Age_binary” (ITS)

View full QIIME2 results (Bray–Curtis – pH_binary)

Figure: PERMANOVA Bray-Curtis for “pH_binary” (ITS)

Differentially Abundant Pathways

To identify pathways with significantly different abundances between groups, ANCOM was applied to PICRUSt2 pathway predictions.

Results Summary

Differential abundance analysis of predicted microbial functions (ANCOM) revealed that several pathways varied significantly in response to environmental factors.

Establishment method strongly influenced functional potential, with two pathways significantly more abundant in specific restoration types: peptidoglycan biosynthesis II (staphylococci) and thiazole component of thiamine diphosphate biosynthesis.

Site age was associated with changes in core metabolic functions, including heme biosynthesis and L-methionine biosynthesis, which were enriched in specific age groups.

Soil pH also shaped functional composition, particularly influencing pathways involved in pentalenolactone and pentalenoketolactone biosynthesis.

Finally, sheep grazing was linked to the greatest number of differentially abundant pathways, notably those involved in glutaryl-CoA degradation, L-glutamate degradation, androstenedione degradation, glycine betaine degradation, and other key metabolic routes.

Select variable: Establishment Cutting Cattle Sheep Plough Age_group pH_binary Age_binary

Differences in Pathway Abundance Associated with 'Establishment' (ANCOM Results) (ITS)

Differences in Pathway Abundance Associated with 'Cutting' (ANCOM Results) (ITS)

Differences in Pathway Abundance Associated with 'Cattle' (ANCOM Results) (ITS)

Differences in Pathway Abundance Associated with 'Sheep' (ANCOM Results) (ITS)

Differences in Pathway Abundance Associated with 'Plough' (ANCOM Results) (ITS)

Differences in Pathway Abundance Associated with 'Age_group' (ANCOM Results) (ITS)

Differences in Pathway Abundance Associated with 'pH_binary' (ANCOM Results) (ITS)

Differences in Pathway Abundance Associated with 'Age_binary' (ANCOM Results) (ITS)

Heatmaps

These heatmaps display the relative abundance of predicted functional pathways across different samples, based on PICRUSt2-inferred metagenomic profiles from 16S rRNA gene data. Generated using QIIME2, each row represents a metabolic pathway (from MetaCyc), and each column corresponds to a sample group. The colour intensity indicates the predicted abundance of each functional pathway in the respective group.

You can click on each image to view a larger version in a new tab.

---

Fungi (ITS)

Results Summary

Fungal functional diversity was primarily influenced by pH and establishment method. Lower functional diversity was observed in alkaline soils. Guild analysis showed dominance of undefined saprotrophs, with dung-associated and pathogenic guilds more prevalent in managed or disturbed sites.

Functional Alpha Diversity

For fungal communities (ITS), functional diversity was assessed using metrics derived from predicted ecological functions (e.g., trophic modes, symbiosis potential).

The metrics below—Shannon diversity and Evenness—summarise the variability of ecological roles across different experimental conditions. Use the dropdown menu to explore patterns by variable.

Results Summary

Establishment method had no significant impact on fungal guild diversity, whether measured by Shannon index or evenness. In contrast, soil pH significantly influenced both metrics, with more acidic to neutral soils (below median pH) showing higher functional diversity than alkaline sites (Shannon: p-adjusted = 0.003; Evenness: p-adjusted = 0.007).

Neither site age nor management practices (e.g., grazing, ploughing, cutting) had a significant effect on functional diversity, although younger sites tended to have slightly lower diversity overall.

Select variable: Establishment Cutting Cattle Sheep Plough Age category pH category Age

Shannon Index Boxplots

Here is a link to the Shannon Index Boxplot for more flexibility on QIIME2: Shannon Index Boxplot - Functional (ITS)

Kruskal-Wallis p-value: 0.529

Kruskal-Wallis p-value: 0.797

Kruskal-Wallis p-value: 0.133

Kruskal-Wallis p-value: 0.592

Kruskal-Wallis p-value: 0.161

Kruskal-Wallis p-value: 0.439

Kruskal-Wallis p-value: 0.00299

Kruskal-Wallis p-value: 0.372

Evenness Index Boxplots

Here is a link to the Evenness Boxplot for more flexibility on QIIME2: Evenness Boxplot - Functional (ITS)

Kruskal-Wallis p-value: 0.226

Kruskal-Wallis p-value: 0.349

Kruskal-Wallis p-value: 0.106

Kruskal-Wallis p-value: 0.617

Kruskal-Wallis p-value: 0.331

Kruskal-Wallis p-value: 0.487

Kruskal-Wallis p-value: 0.00668

Kruskal-Wallis p-value: 0.206

Functional Beta Diversity

Results Summary

Beta diversity analysis of fungal guild composition revealed that establishment method, soil pH, and site age significantly influenced community structure, as assessed by both Bray–Curtis and Jaccard distance metrics.

PCoA plots showed some clustering by establishment type. Sites using green hay or seed mix were more widely dispersed, while natural regeneration (NR) sites formed a tighter and more cohesive cluster. PERMANOVA confirmed significant differences between establishment types for both Bray–Curtis (p = 0.003) and Jaccard (p = 0.002), with all pairwise comparisons significant. The strongest contrasts were observed between GH and NR for Bray–Curtis and between NR vs. GH/SM for Jaccard.

For pH category, clustering was less distinct with Bray–Curtis but more apparent with Jaccard, where sites with slightly acidic to neutral soils grouped more closely. PERMANOVA results supported these observations, with stronger significance for Jaccard (p-adjusted = 0.001) than for Bray–Curtis (p-adjusted = 0.003).

Site age also had a significant effect on guild composition for both metrics. The effect was more pronounced with Bray–Curtis (p-adjusted = 0.009) than with Jaccard (p-adjusted = 0.04). While PCoA plots suggested some separation by age, each age group contained a mix of sites, indicating partial overlap in guild structure across age categories.

Bray-Curtis

Emperor Plot

Here is a link to the Bray-Curtis Emperor Plot for more flexibility on QIIME2: Bray-Curtis Emperor Plot (16S)

Bray-Curtis Emperor Plot (Fungi - ITS, Functional)

PERMANOVA

You can click on the images below to access the full QIIME2 report.

Select variable: Establishment pH category Age category

View full QIIME2 results (Bray–Curtis – Establishment)

Figure: PERMANOVA Bray-Curtis for “Establishment” (ITS)

View full QIIME2 results (Bray–Curtis – pH)

Figure: PERMANOVA Bray-Curtis for “pH” (ITS)

View full QIIME2 results (Bray–Curtis – Age)

Figure: PERMANOVA Bray-Curtis for “Age” (ITS)

Jaccard

Emperor Plot

Here is a link to the Jaccard Emperor Plot for more flexibility on QIIME2: Bray-Curtis Emperor Plot (16S)

Jaccard Emperor Plot (Fungi - ITS, Functional)

PERMANOVA

You can click on the images below to access the full QIIME2 report.

Select variable: Establishment pH category Age category

View full QIIME2 results (Jaccard – Establishment)

Figure: PERMANOVA Jaccard for “Establishment” (ITS)

View full QIIME2 results (Jaccard – pH)

Figure: PERMANOVA Jaccard for “pH” (ITS)

View full QIIME2 results (Jaccard – Age)

Figure: PERMANOVA Jaccard for “Age” (ITS)

Fungal Functional Guilds (FUNGuild)

In microbial ecology, a guild refers to a group of organisms that fulfil similar ecological roles, regardless of their taxonomic identity. Understanding functional guilds allows researchers to move beyond taxonomic profiles and assess the ecological roles that microbial communities may play in an environment.

To investigate the ecological roles of fungal communities, we used FUNGuild, a tool that assigns fungi to functional guilds based on curated databases and literature. These guilds represent ecological strategies such as:

• Saprotrophs: decomposers of organic matter
  
• Mycorrhizal fungi: symbionts associated with plant roots
  
• Pathogens: organisms that cause disease in plants or animals
• Ectomycorrhizal: symbionts that assist plant roots by enhancing nutrient and water uptake, especially nitrogen and phosphorus, while contributing to carbon cycling and soil nutrient mobilisation
• Arbuscular Mycorrhizal: symbionts that contribute to nutrient uptake
• Endophyte: microorganisms that promote the growth and development of the plants
• Lichenized: fungi that form symbiotic associations with algae or cyanobacteria, creating lichens that can survive in harsh environments by combining structural support and photosynthetic ability
• Parasites: fungi that live on or inside a host organism, extracting nutrients and often causing harm
• Symbiotrophs: fungi that engage in mutually beneficial relationships with host organisms, typically exchanging nutrients for resources like carbon or protection

This functional classification provides valuable insights into what fungi are likely doing in the ecosystem, beyond simply who they are.

This section explores the functional roles of fungi within each site, based on guild-level annotations provided by FUNGuild. Fungal guilds reflect ecological functions such as saprotrophy, symbiosis (e.g., mycorrhizal fungi), or pathogenicity. This approach provides insight into how fungal communities may contribute to ecosystem processes, complementing traditional taxonomic analyses.

Top 20 Most Abundant Guilds

The plot below highlights the top 20 most abundant fungal guilds identified using FUNGuild. To avoid clutter, the guild names are hidden on the y-axis; however, users can hover over each bar to reveal the full name, enabling interactive and detailed exploration of fungal functional diversity.

Top 20 functional guilds

Guild Abundance Across Sites

The figure below shows the total abundance of fungal ASVs across sites, aggregated by functional guild. This provides an overview of how guild-level composition varies between locations, which may reflect differences in land use, soil conditions, or restoration histories.

Fungal communities were dominated by Undefined Saprotrophs, Plant Pathogens, and Dung Saprotrophs. A high proportion of Unassigned fungi highlights limits in current taxonomic resolution.

Fungal Guild Abundance by Site

Guild Abundance by Variable

Establishment

Natural regeneration had the highest proportion of Undefined Saprotrophs, while Green Hay/Bush and Seed Mix had more Plant Pathogens and Dung Saprotrophs.

pH

Alkaline sites showed dominance of Undefined Saprotrophs; acidic sites had more Animal Parasites.

Age

Older sites had more Undefined Saprotrophs; younger sites showed higher proportions of Plant Pathogens and Dung Saprotrophs.

Select variable: Establishment Age pH

Fungal Guild Relative Abundance by pH_category

Fungal Guild Relative Abundance by Age_category

Fungal Guild Relative Abundance by Establishment

---

Bacterial and Fungal Comparison

Results Summary
Correlations between bacterial and fungal taxa and functions revealed potential ecological links, particularly between saprotrophic/pathogenic fungi and sugar-degrading or nitrogen-associated bacterial pathways. Most associations were site-specific and stronger in NR plots, with nutrient cycling emerging as a shared functional axis.

Taxonomic level

After filtering for the most abundant taxa, seven significant correlations were detected at the family level, with Pearson correlation coefficients ranging from -0.61 to 0.80. Among these, six were positive and one was negative.

The main bacterial families involved were:

• Bradyrhizobiaceae
• Thermomonosporaceae
• Rhodospirillaceae
• Chthoniobacteraceae
• Rhodobacteraceae

The corresponding fungal families included:

• Trimorphomycetaceae
• Clavicipitaceae
• Tricholomataceae
• Piskurozymaceae
• Pseudeurotiaceae

At the class level, dominant bacterial groups were Alphaproteobacteria, Actinobacteria, and Spartobacteria. Fungal communities were mainly composed of Tremellomycetes, Sordariomycetes, Agaricomycetes, and Leotiomycetes.

Below is an interactive visualisation of significant correlations between abundant bacterial and fungal families:

Click to view the detailed correlation table

Significant correlations between bacterial and fungal taxa

Bacteria	Fungi	r	p_adj
Bradyrhizobiaceae	Trimorphomycetaceae	0.80	0.0000000
Bradyrhizobiaceae	Clavicipitaceae	0.79	0.0000000
Thermomonosporaceae	Trimorphomycetaceae	0.73	0.0000007
Rhodospirillaceae	Tricholomataceae	0.72	0.0000022
Thermomonosporaceae	Piskurozymaceae	0.72	0.0000028
Chthoniobacteraceae	Pseudeurotiaceae	0.72	0.0000039
Rhodobacteraceae	Clavicipitaceae	-0.61	0.0118590

---

Functional level

A total of 42 significant correlations were identified between bacterial pathways (from the MetaCyc database) and fungal guilds (from FUNGuild). Among these, 40 correlations were positive, and 2 were negative, with Pearson correlation values ranging from -0.75 to 0.98.

While most correlations involved rare pathways or guilds (appearing at only one or two sites or with low overall counts), 10 correlations stood out due to their association with more abundant functions.

The dominant fungal guilds among these included:

• Animal pathogen / dung saprotrophs
• Endophytic / ericoid saprotrophs
• Fungal and insect parasites

The correlated bacterial pathways primarily involved sugar and protein degradation, such as:

• Methylglyoxal degradation
• L-arabinose degradation
• Glucose, sucrose, and L-tryptophan degradation
• Purine degradation

Below is an interactive visualisation of significant correlations between abundant bacterial pathways and fungal guilds:

Click to view the detailed table for the 10 strongest correlations

Significant correlations between bacterial pathways and fungal guilds

Guild	Pathway	r	p_ajd
Animal pathogen / dung saprotroph	Methylglyoxal degradation	0.94	0.0e+00
Animal pathogen / dung saprotroph	L-arabinose degradation IV	0.85	0.0e+00
Animal pathogen / dung saprotroph	Glucose and glucose-1-phosphate degradation	0.80	0.0e+00
Animal pathogen / dung saprotroph	Sucrose degradation IV (sucrose phosphorylase)	0.78	0.0e+00
Endophyte–Ericoid Mycorrhizal–Undefined Saprotroph	Superpathway of N-acetylneuraminate, N-acetylglucosamine, and N-acetylmannosamine degradation to β-D-fructofuranose 6-phosphate	0.74	1.0e-07
Animal pathogen / dung saprotroph	L-tryptophan degradation to 2-amino-3-carboxymuconate semialdehyde	0.72	8.0e-07
Animal pathogen / dung saprotroph	2-methylcitrate cycle II	0.71	1.4e-06
Endophyte–Ericoid Mycorrhizal–Undefined Saprotroph	Superpathway of N-acetylneuraminate degradation	0.71	2.4e-06
Dung Saprotroph–Ectomycorrhizal–Litter Saprotroph	Purine nucleobases degradation I (anaerobic)	-0.71	2.6e-06
Endophyte–Epiphyte–Fungal Parasite–Insect Parasite	Nylon-6 oligomer degradation	0.71	3.0e-06

Discussion

Bacterial Alpha and Beta Diversity

Microbial alpha diversity was significantly higher in assisted restoration plots (GH and SM) compared to natural regeneration (NR), suggesting that sown vegetation enhances microbial richness and evenness, likely via increased plant inputs and rhizosphere complexity (Van der Putten et al., 2013; Bardgett & van der Wal, 2014). In contrast, soil pH, management, and site age had no significant impact on alpha diversity, possibly due to functional redundancy that stabilises core microbial functions (Louca et al., 2018). Beta diversity was strongly shaped by all tested factors—restoration method, site age, soil pH, and grazing—particularly establishment method, which significantly altered both community membership (Unweighted UniFrac, p = 0.001) and abundance profiles (Bray-Curtis, p = 0.001). Site age influenced Bray-Curtis dissimilarities but not UniFrac distances, suggesting abundance shifts without major taxonomic turnover (Fierer et al., 2010; Nemergut et al., 2013). Soil pH acted as a strong environmental filter, affecting microbial community structure by influencing enzyme activity, nutrient availability, and stress response (Lauber et al., 2009; Rousk et al., 2010). Grazing also altered beta diversity, likely through disturbance and nutrient inputs (Bardgett & Wardle, 2010; Patra et al., 2005).

Functional Diversity and Key Drivers

Functional beta diversity was also significantly affected by environmental gradients, with the restoration method again exerting the strongest influence (PERMANOVA p = 0.001). Functional divergence likely reflected changes in rhizosphere chemistry and substrate input driven by vegetation type (Eisenhauer et al., 2017; Zhalnina et al., 2018; Lange et al., 2015). Soil pH had a strong effect on functional profiles (PERMANOVA p = 0.012), with alkaline soils enriched in pathways related to secondary metabolism (e.g., pentalenolactone biosynthesis), a trait associated with antibiotic production in Actinobacteria (Tetzlaff et al., 2006). Site age had a moderate effect (p = 0.035), with younger sites enriched in biosynthetic pathways such as heme and methionine biosynthesis (Choby & Skaar, 2016; Mohany et al., 2021). These functions are linked to early successional microbial activity (Fierer et al., 2010; Nemergut et al., 2013). Grazing influenced functions tied to nitrogen and organic matter degradation (Wang et al., 2020; Tang et al., 2020).

Differential Pathways and Functional Responses

Pathway-level differences revealed microbial adaptation to specific environmental contexts. Alkaline soils were enriched in antibiotic-related pathways (Tetzlaff et al., 2006), while younger sites showed prioritisation of growth-related functions (Choby & Skaar, 2016). In GH and SM plots, increased peptidoglycan and thiamine biosynthesis pathways suggested higher microbial growth and activity (Heijenoort, 2001; Begley et al., 1999; Eisenhauer et al., 2017). These adaptations reflect changes in vegetation and nutrient input linked to restoration strategy.

Core Metabolism and Redundancy

Despite compositional variation, core biosynthetic and energy metabolism pathways were conserved across all sites, including amino acid synthesis and pyruvate fermentation. This points to the persistence of essential ecosystem functions and supports the concept of functional redundancy (Louca et al., 2018).

Taxonomic Trends and Indicators

Dominant phyla included Actinobacteria and Proteobacteria, consistent with nutrient-rich conditions in SM and GH plots (Spain et al., 2009). Acidobacteria, associated with acidic soils, remained common even in slightly acidic environments (Jones et al., 2009). Conexibacteraceae were more abundant in lower-pH soils (Lauber et al., 2009; Fierer & Jackson, 2006), while Intrasporangiaceae and Blastococcus preferred SM and GH sites, responding to increased organic inputs (Jangid et al., 2011).

Fungal Diversity and Environmental Drivers

Fungal diversity was mainly shaped by pH and site age. Alpha diversity decreased in alkaline soils, potentially due to reduced nutrient use efficiency in fungi (Ma et al., 2023), contradicting previous findings of increased diversity at high pH (Seaton et al., 2022). Older sites supported more complex fungal communities, consistent with longer recovery times post-disturbance (Gao et al., 2021). Establishment method also influenced fungal composition, with GH and SM sites hosting distinct taxa.

Fungal Composition and Functional Guilds

Metarhizium was more frequent in NR sites, reflecting reduced disturbance (Zimmermann, 2007). Gibellulopsis, potentially linked to pathogenic interactions, appeared in GH and SM plots (Alisaac & Götz, 2022). Sheep grazing was associated with higher Helotiales and dung-associated fungi like Preussia (Kruys et al., 2006; Bruyant et al., 2024), while ploughing promoted taxa like Saitozyma capable of breaking down plant-derived carbon (Aliyu et al., 2021). Site age influenced composition, with early successional fungi like Lasiosphaeriaceae more common in younger plots, and potential endophytes like Paraphaeosphaeria found in older sites (Baroncelli et al., 2020; Louw et al., 2022).

Fungal Functional Guild Shifts

Functional guild analysis showed lower abundance of pathogens and dung saprotrophs in NR sites, with higher prevalence in GH and SM plots. pH subtly influenced parasite abundance, while site age showed more complex guild structuring, aligning with patterns observed in forest restoration studies (Fang et al., 2023). Across all sites, saprotrophs dominated, particularly unclassified ones, underlining gaps in fungal ecological data (Nilsson et al., 2020).

Fungal-Bacterial Comparison and Interactions While bacterial diversity was most strongly shaped by establishment method, fungal diversity was more affected by pH and age. However, community composition (beta diversity) for both kingdoms was driven by similar environmental filters. Cross-kingdom correlations included co-occurrence of taxa in GH and SM sites, such as Proteobacteria and Gibellulopsis. Functional comparisons showed complementarity in decomposition and nutrient cycling, supporting previous observations of cross-kingdom cooperation (Wagg et al., 2019; Kohlmeier et al., 2005).

Functional and Taxonomic Correlations

Bradyrhizobiaceae (oligotrophic, nitrogen-fixing bacteria) were associated with saprotrophic yeasts (Saitozyma) and parasitic fungi (Keithomyces), possibly due to shared preference for low pH (Sawada et al., 2023). Rhodospirillaceae and Dermoloma co-occurred on one site, suggesting shared habitat preferences. Most functional correlations involved nutrient cycling: Thermothielavioides (a dung saprotroph) was associated with sugar degradation pathways and propionate metabolism (Camargo et al., 2020; Caspi et al., 2020). These patterns indicate that bacteria and fungi carry out complementary metabolic roles, particularly in NR sites, where successional complexity may support tighter functional coupling (Gao et al., 2021).

Limitations

Correlation analyses were limited by the complexity and abundance of fungal-bacterial interactions. Filtering was required to manage the volume of pairwise comparisons, which may have excluded meaningful associations. The ambiguity of FUNGuild annotations also constrained ecological interpretation. Alternative modelling approaches like GLMs may yield more robust insights into cross-kingdom interactions. Additionally, the large proportion of unassigned fungi also underscores the limitations of current fungal functional classification systems (Nilsson et al., 2020), potentially obscuring a portion of the ecological signal.